library(tidyverse)
library(plotly)
library(sf)
library(mapview)
library(tigris)
library(censusapi)
library(leaflet)
library(lehdr)
library(usmap)
library(lmtest)
library(pracma)
library(lmtest)
library(forecast)
library(vars)
library(rvest)
library(RSelenium)
library(seleniumPipes)
library(dLagM)
options(
tigris_class = "sf",
tigris_use_cache = TRUE
)
Get the case data, first for SCC.
# remDr <- remoteDriver(
# remoteServerAddr = "192.168.86.25",
# port = 4445L
# )
# remDr$open()
#
# remDr$navigate("https://app.powerbigov.us/view?r=eyJrIjoiZTg2MTlhMWQtZWE5OC00ZDI3LWE4NjAtMTU3YWYwZDRlOTNmIiwidCI6IjBhYzMyMDJmLWMzZTktNGY1Ni04MzBkLTAxN2QwOWQxNmIzZiJ9")
#
# webElem <- remDr$findElements(using = "class", value = "column")
#
# cases <-
# 1:length(webElem) %>%
# map(function(x){
# webElem[[x]]$getElementAttribute("aria-label") %>% as.character()
# }) %>%
# unlist() %>%
# as.data.frame()
#
# scc_cumulative_cases <-
# cases %>%
# rename(text = ".") %>%
# filter(grepl("Total_cases",text)) %>%
# separate(text, c("date","cases"), sep = "\\.") %>%
# mutate(
# date =
# substr(date,6,nchar(.)) %>%
# as.Date("%A, %B %d, %Y"),
# cases =
# substr(cases,13,nchar(.)) %>%
# as.numeric()
# )
#
# saveRDS(scc_cumulative_cases, file = "/Users/simonespeizer/Documents/2020 Spring Quarter/CEE 218Z/covid19/Simone_Speizer/scc_cumulative_cases.rds")
scc_cumulative_cases <- readRDS("/Users/simonespeizer/Documents/2020 Spring Quarter/CEE 218Z/covid19/Simone_Speizer/scc_cumulative_cases.rds")
Also for SMC.
# remDr$navigate("https://app.powerbigov.us/view?r=eyJrIjoiMWI5NmE5M2ItOTUwMC00NGNmLWEzY2UtOTQyODA1YjQ1NWNlIiwidCI6IjBkZmFmNjM1LWEwNGQtNDhjYy1hN2UzLTZkYTFhZjA4ODNmOSJ9")
#
# webElem <- remDr$findElements(using = "class", value = "column")
#
# tests <-
# 1:length(webElem) %>%
# map(function(x){
# webElem[[x]]$getElementAttribute("aria-label") %>% as.character()
# }) %>%
# unlist() %>%
# as.data.frame()
#
# tests_clean <-
# tests %>%
# rename(text = ".") %>%
# separate(text, c("date","test_text"), sep = "\\.") %>%
# separate(test_text, c(NA,"type",NA,"tests")) %>%
# mutate(
# date =
# substr(date,23,nchar(.)) %>%
# as.Date("%A, %B %d, %Y"),
# tests =
# tests %>%
# as.numeric()
# ) %>%
# spread(
# key = type,
# value = tests
# ) %>%
# mutate(
# total = Negative + Pending + Positive,
# perc_positive = Positive/total,
# perc_positive_movavg = movavg(perc_positive, 7, type = "s")
# )
#
# smc_cumulative_cases <- tests_clean %>%
# mutate(cumulative_cases = cumsum(Positive), cumulative_negative = cumsum(Negative), cumulative_total = cumulative_cases+cumulative_negative, perc_positive_cumulative = cumulative_cases*100 / cumulative_total, perc_positive_cumulative_mov = movavg(perc_positive_cumulative, 7, type = "s"))
#
# saveRDS(smc_cumulative_cases, file = "/Users/simonespeizer/Documents/2020 Spring Quarter/CEE 218Z/covid19/Simone_Speizer/smc_cumulative_cases.rds")
smc_cumulative_cases <- readRDS("/Users/simonespeizer/Documents/2020 Spring Quarter/CEE 218Z/covid19/Simone_Speizer/smc_cumulative_cases.rds")
Get social distancing data.
scc_blockgroups <-
block_groups("CA","Santa Clara", cb=F, progress_bar=F) %>%
st_transform('+proj=longlat +datum=WGS84')
smc_blockgroups <-
block_groups("CA","San Mateo", cb=F, progress_bar=F) %>%
st_transform('+proj=longlat +datum=WGS84')
bay_sd <- readRDS("/Users/simonespeizer/pCloud Drive/Shared/SFBI/Restricted Data Library/Safegraph/covid19analysis/bay_socialdistancing_v2.rds") %>%
mutate(date = date_range_start %>% substr(1,10) %>% as.Date())
# obtaining weekends
weekends <- bay_sd %>%
filter(!duplicated(date)) %>%
arrange(date) %>%
mutate(weekend = ifelse((date %>% as.numeric()) %% 7 %in% c(2,3), T, F)) %>%
dplyr::select(date,weekend)
bay_sd <- bay_sd %>% left_join(weekends)
scc_cases_sd_daily <- bay_sd %>%
filter(origin_census_block_group %in% scc_blockgroups$GEOID) %>%
group_by(date) %>%
summarize(total_at_home = sum(completely_home_device_count), total_devices = sum(device_count)) %>%
mutate(
percent_at_home = total_at_home*100/total_devices,
percent_leaving_home = (100 - percent_at_home),
) %>%
left_join(
scc_cumulative_cases
) %>%
filter(date >= min(scc_cumulative_cases$date))
# get the baseline percent of people at home
pre_case_growth_at_home_scc <- bay_sd %>%
filter(date < min(scc_cumulative_cases$date)) %>%
filter(origin_census_block_group %in% scc_blockgroups$GEOID) %>%
summarize(total_at_home = sum(completely_home_device_count), total_devices = sum(device_count)) %>%
mutate(percent_at_home = total_at_home*100/total_devices, percent_leaving_home = (100 - percent_at_home))
# include change in percent leaving home
scc_cases_sd_daily <- scc_cases_sd_daily %>%
mutate(leaving_home_dif = percent_leaving_home - pre_case_growth_at_home_scc$percent_leaving_home[1],
log_cases = log(cases))
# compute number of differences for stationarity
ndiffs(scc_cases_sd_daily$cases)
## [1] 2
ndiffs(scc_cases_sd_daily$log_cases[-1])
## [1] 2
ndiffs(scc_cases_sd_daily$leaving_home_dif)
## [1] 1
scc_test_big <- scc_cases_sd_daily %>%
mutate(
dlog_cases = c(NA,diff(log_cases)),
d2log_cases = c(NA,NA,diff(log_cases, differences = 2)),
dcases = c(NA,diff(cases)),
d2cases = c(NA,NA,diff(dcases, differences = 2)),
dleaving = c(NA,diff(leaving_home_dif)),
d2leaving = c(NA,NA,diff(leaving_home_dif, differences = 2)),
cases_mov = movavg(cases, 7, type = "s"),
log_cases_mov = movavg(log_cases, 7, type = "s"),
dlog_cases_mov = c(NA,diff(log_cases_mov)),
d2log_cases_mov = c(NA,NA,diff(log_cases_mov, differences = 2)),
dcases_mov = c(NA,diff(cases_mov)),
d2cases_mov = c(NA,diff(dcases_mov)),
leaving_mov = movavg(leaving_home_dif, 7, type = "s"),
dleaving_mov = c(NA,diff(leaving_mov)),
d2leaving_mov = c(NA,diff(dleaving_mov)),
leaving4 = c(rep(NA,4), scc_cases_sd_daily$leaving_home_dif[1:(nrow(scc_cases_sd_daily)-4)]),
dleaving4 = c(NA,diff(leaving4)),
d2leaving4 = c(NA,NA,diff(leaving4, differences = 2)),
leaving3 = c(rep(NA,3), scc_cases_sd_daily$leaving_home_dif[1:(nrow(scc_cases_sd_daily)-3)]),
leaving3_mov = movavg(leaving3, 7, type = "s"),
dleaving3_mov = c(NA,diff(leaving3_mov)),
d2leaving3_mov = c(NA,NA,diff(leaving3_mov, differences = 2)),
leaving4_mov = movavg(leaving4, 7, type = "s"),
dleaving4_mov = c(NA,diff(leaving4_mov)),
d2leaving4_mov = c(NA,NA,diff(leaving4_mov, differences = 2)),
leaving5 = c(rep(NA,5), scc_cases_sd_daily$leaving_home_dif[1:(nrow(scc_cases_sd_daily)-5)]),
leaving5_mov = movavg(leaving5, 7, type = "s"),
dleaving5_mov = c(NA,diff(leaving5_mov)),
d2leaving5_mov = c(NA,NA,diff(leaving5_mov, differences = 2)),
leaving6 = c(rep(NA,6), scc_cases_sd_daily$leaving_home_dif[1:(nrow(scc_cases_sd_daily)-6)]),
leaving6_mov = movavg(leaving6, 7, type = "s"),
dleaving6_mov = c(NA,diff(leaving6_mov)),
d2leaving6_mov = c(NA,NA,diff(leaving6_mov, differences = 2)),
leaving7 = c(rep(NA,7), scc_cases_sd_daily$leaving_home_dif[1:(nrow(scc_cases_sd_daily)-7)]),
leaving7_mov = movavg(leaving7, 7, type = "s"),
dleaving7_mov = c(NA,diff(leaving7_mov)),
d2leaving7_mov = c(NA,NA,diff(leaving7_mov, differences = 2)),
leaving8 = c(rep(NA,8), scc_cases_sd_daily$leaving_home_dif[1:(nrow(scc_cases_sd_daily)-8)]),
leaving8_mov = movavg(leaving8, 7, type = "s"),
dleaving8_mov = c(NA,diff(leaving8_mov)),
d2leaving8_mov = c(NA,NA,diff(leaving8_mov, differences = 2)),
leaving9 = c(rep(NA,9), scc_cases_sd_daily$leaving_home_dif[1:(nrow(scc_cases_sd_daily)-9)]),
leaving9_mov = movavg(leaving9, 7, type = "s"),
dleaving9_mov = c(NA,diff(leaving9_mov)),
d2leaving9_mov = c(NA,NA,diff(leaving9_mov, differences = 2)),
leaving10 = c(rep(NA,10), scc_cases_sd_daily$leaving_home_dif[1:(nrow(scc_cases_sd_daily)-10)]),
leaving10_mov = movavg(leaving10, 7, type = "s"),
dleaving10_mov = c(NA,diff(leaving10_mov)),
d2leaving10_mov = c(NA,NA,diff(leaving10_mov, differences = 2)),
past_cases = c(NA, scc_cases_sd_daily$cases[1:(nrow(scc_cases_sd_daily)-1)]),
cases_growth_daily = (dcases / past_cases)*100,
cases_growth_daily_mov = movavg(cases_growth_daily, 7, type = "s")
) %>%
filter(date >= "2020-03-01")
ndiffs(scc_cases_sd_daily$cases_growth_daily)
## [1] 0
Granger
grangertest(
d2log_cases_mov ~ dleaving_mov,
order = 3,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:3) + Lags(dleaving_mov, 1:3)
## Model 2: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:3)
## Res.Df Df F Pr(>F)
## 1 73
## 2 76 -3 9.0988 3.424e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving_mov ~ d2log_cases_mov,
order = 3,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving_mov ~ Lags(dleaving_mov, 1:3) + Lags(d2log_cases_mov, 1:3)
## Model 2: dleaving_mov ~ Lags(dleaving_mov, 1:3)
## Res.Df Df F Pr(>F)
## 1 73
## 2 76 -3 1.2405 0.3013
grangertest(
d2log_cases_mov ~ dleaving_mov,
order = 4,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:4) + Lags(dleaving_mov, 1:4)
## Model 2: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:4)
## Res.Df Df F Pr(>F)
## 1 70
## 2 74 -4 8.1953 1.764e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving_mov ~ d2log_cases_mov,
order = 4,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving_mov ~ Lags(dleaving_mov, 1:4) + Lags(d2log_cases_mov, 1:4)
## Model 2: dleaving_mov ~ Lags(dleaving_mov, 1:4)
## Res.Df Df F Pr(>F)
## 1 70
## 2 74 -4 0.8554 0.4951
grangertest(
d2log_cases_mov ~ dleaving_mov,
order = 5,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:5) + Lags(dleaving_mov, 1:5)
## Model 2: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:5)
## Res.Df Df F Pr(>F)
## 1 67
## 2 72 -5 6.7985 3.506e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving_mov ~ d2log_cases_mov,
order = 5,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving_mov ~ Lags(dleaving_mov, 1:5) + Lags(d2log_cases_mov, 1:5)
## Model 2: dleaving_mov ~ Lags(dleaving_mov, 1:5)
## Res.Df Df F Pr(>F)
## 1 67
## 2 72 -5 0.694 0.6298
grangertest(
d2log_cases_mov ~ dleaving_mov,
order = 6,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:6) + Lags(dleaving_mov, 1:6)
## Model 2: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:6)
## Res.Df Df F Pr(>F)
## 1 64
## 2 70 -6 5.3088 0.0001703 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving_mov ~ d2log_cases_mov,
order = 6,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving_mov ~ Lags(dleaving_mov, 1:6) + Lags(d2log_cases_mov, 1:6)
## Model 2: dleaving_mov ~ Lags(dleaving_mov, 1:6)
## Res.Df Df F Pr(>F)
## 1 64
## 2 70 -6 1.041 0.4074
grangertest(
d2log_cases_mov ~ dleaving_mov,
order = 7,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:7) + Lags(dleaving_mov, 1:7)
## Model 2: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:7)
## Res.Df Df F Pr(>F)
## 1 61
## 2 68 -7 2.9586 0.009776 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving_mov ~ d2log_cases_mov,
order = 7,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving_mov ~ Lags(dleaving_mov, 1:7) + Lags(d2log_cases_mov, 1:7)
## Model 2: dleaving_mov ~ Lags(dleaving_mov, 1:7)
## Res.Df Df F Pr(>F)
## 1 61
## 2 68 -7 1.8237 0.09882 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
d2log_cases_mov ~ dleaving_mov,
order = 8,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:8) + Lags(dleaving_mov, 1:8)
## Model 2: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:8)
## Res.Df Df F Pr(>F)
## 1 58
## 2 66 -8 4.8424 0.0001338 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving_mov ~ d2log_cases_mov,
order = 8,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving_mov ~ Lags(dleaving_mov, 1:8) + Lags(d2log_cases_mov, 1:8)
## Model 2: dleaving_mov ~ Lags(dleaving_mov, 1:8)
## Res.Df Df F Pr(>F)
## 1 58
## 2 66 -8 2.1953 0.04087 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
d2log_cases_mov ~ dleaving_mov,
order = 9,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:9) + Lags(dleaving_mov, 1:9)
## Model 2: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:9)
## Res.Df Df F Pr(>F)
## 1 55
## 2 64 -9 3.851 0.0007892 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving_mov ~ d2log_cases_mov,
order = 9,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving_mov ~ Lags(dleaving_mov, 1:9) + Lags(d2log_cases_mov, 1:9)
## Model 2: dleaving_mov ~ Lags(dleaving_mov, 1:9)
## Res.Df Df F Pr(>F)
## 1 55
## 2 64 -9 1.7838 0.09234 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
d2log_cases_mov ~ dleaving_mov,
order = 10,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:10) + Lags(dleaving_mov, 1:10)
## Model 2: d2log_cases_mov ~ Lags(d2log_cases_mov, 1:10)
## Res.Df Df F Pr(>F)
## 1 52
## 2 62 -10 4.0521 0.0003831 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving_mov ~ d2log_cases_mov,
order = 10,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving_mov ~ Lags(dleaving_mov, 1:10) + Lags(d2log_cases_mov, 1:10)
## Model 2: dleaving_mov ~ Lags(dleaving_mov, 1:10)
## Res.Df Df F Pr(>F)
## 1 52
## 2 62 -10 1.5036 0.1648
Granger for not moving average
grangertest(
d2log_cases ~ dleaving,
order = 3,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases ~ Lags(d2log_cases, 1:3) + Lags(dleaving, 1:3)
## Model 2: d2log_cases ~ Lags(d2log_cases, 1:3)
## Res.Df Df F Pr(>F)
## 1 73
## 2 76 -3 4.467 0.006171 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ d2log_cases,
order = 3,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:3) + Lags(d2log_cases, 1:3)
## Model 2: dleaving ~ Lags(dleaving, 1:3)
## Res.Df Df F Pr(>F)
## 1 73
## 2 76 -3 0.9464 0.4228
grangertest(
d2log_cases ~ dleaving,
order = 4,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases ~ Lags(d2log_cases, 1:4) + Lags(dleaving, 1:4)
## Model 2: d2log_cases ~ Lags(d2log_cases, 1:4)
## Res.Df Df F Pr(>F)
## 1 70
## 2 74 -4 6.5116 0.0001631 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ d2log_cases,
order = 4,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:4) + Lags(d2log_cases, 1:4)
## Model 2: dleaving ~ Lags(dleaving, 1:4)
## Res.Df Df F Pr(>F)
## 1 70
## 2 74 -4 1.24 0.302
grangertest(
d2log_cases ~ dleaving,
order = 5,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases ~ Lags(d2log_cases, 1:5) + Lags(dleaving, 1:5)
## Model 2: d2log_cases ~ Lags(d2log_cases, 1:5)
## Res.Df Df F Pr(>F)
## 1 67
## 2 72 -5 4.1674 0.002344 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ d2log_cases,
order = 5,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:5) + Lags(d2log_cases, 1:5)
## Model 2: dleaving ~ Lags(dleaving, 1:5)
## Res.Df Df F Pr(>F)
## 1 67
## 2 72 -5 0.8606 0.5123
grangertest(
d2log_cases ~ dleaving,
order = 6,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases ~ Lags(d2log_cases, 1:6) + Lags(dleaving, 1:6)
## Model 2: d2log_cases ~ Lags(d2log_cases, 1:6)
## Res.Df Df F Pr(>F)
## 1 64
## 2 70 -6 2.9787 0.01253 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ d2log_cases,
order = 6,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:6) + Lags(d2log_cases, 1:6)
## Model 2: dleaving ~ Lags(dleaving, 1:6)
## Res.Df Df F Pr(>F)
## 1 64
## 2 70 -6 0.844 0.5408
grangertest(
d2log_cases ~ dleaving,
order = 7,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases ~ Lags(d2log_cases, 1:7) + Lags(dleaving, 1:7)
## Model 2: d2log_cases ~ Lags(d2log_cases, 1:7)
## Res.Df Df F Pr(>F)
## 1 61
## 2 68 -7 2.7765 0.01423 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ d2log_cases,
order = 7,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:7) + Lags(d2log_cases, 1:7)
## Model 2: dleaving ~ Lags(dleaving, 1:7)
## Res.Df Df F Pr(>F)
## 1 61
## 2 68 -7 0.9196 0.4978
grangertest(
d2log_cases ~ dleaving,
order = 8,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases ~ Lags(d2log_cases, 1:8) + Lags(dleaving, 1:8)
## Model 2: d2log_cases ~ Lags(d2log_cases, 1:8)
## Res.Df Df F Pr(>F)
## 1 58
## 2 66 -8 2.765 0.01167 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ d2log_cases,
order = 8,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:8) + Lags(d2log_cases, 1:8)
## Model 2: dleaving ~ Lags(dleaving, 1:8)
## Res.Df Df F Pr(>F)
## 1 58
## 2 66 -8 1.2574 0.2835
grangertest(
d2log_cases ~ dleaving,
order = 9,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases ~ Lags(d2log_cases, 1:9) + Lags(dleaving, 1:9)
## Model 2: d2log_cases ~ Lags(d2log_cases, 1:9)
## Res.Df Df F Pr(>F)
## 1 55
## 2 64 -9 2.712 0.01089 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ d2log_cases,
order = 9,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:9) + Lags(d2log_cases, 1:9)
## Model 2: dleaving ~ Lags(dleaving, 1:9)
## Res.Df Df F Pr(>F)
## 1 55
## 2 64 -9 1.4197 0.2025
grangertest(
d2log_cases ~ dleaving,
order = 10,
data = scc_test_big
)
## Granger causality test
##
## Model 1: d2log_cases ~ Lags(d2log_cases, 1:10) + Lags(dleaving, 1:10)
## Model 2: d2log_cases ~ Lags(d2log_cases, 1:10)
## Res.Df Df F Pr(>F)
## 1 52
## 2 62 -10 2.2658 0.02756 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ d2log_cases,
order = 10,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:10) + Lags(d2log_cases, 1:10)
## Model 2: dleaving ~ Lags(dleaving, 1:10)
## Res.Df Df F Pr(>F)
## 1 52
## 2 62 -10 0.9806 0.4715
Granger with case growth rate
grangertest(
cases_growth_daily ~ dleaving,
order = 2,
data = scc_test_big
)
## Granger causality test
##
## Model 1: cases_growth_daily ~ Lags(cases_growth_daily, 1:2) + Lags(dleaving, 1:2)
## Model 2: cases_growth_daily ~ Lags(cases_growth_daily, 1:2)
## Res.Df Df F Pr(>F)
## 1 76
## 2 78 -2 0.8033 0.4516
grangertest(
dleaving ~ cases_growth_daily,
order = 2,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:2) + Lags(cases_growth_daily, 1:2)
## Model 2: dleaving ~ Lags(dleaving, 1:2)
## Res.Df Df F Pr(>F)
## 1 76
## 2 78 -2 3.6329 0.03113 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
cases_growth_daily ~ dleaving,
order = 3,
data = scc_test_big
)
## Granger causality test
##
## Model 1: cases_growth_daily ~ Lags(cases_growth_daily, 1:3) + Lags(dleaving, 1:3)
## Model 2: cases_growth_daily ~ Lags(cases_growth_daily, 1:3)
## Res.Df Df F Pr(>F)
## 1 73
## 2 76 -3 2.2196 0.09306 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ cases_growth_daily,
order = 3,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:3) + Lags(cases_growth_daily, 1:3)
## Model 2: dleaving ~ Lags(dleaving, 1:3)
## Res.Df Df F Pr(>F)
## 1 73
## 2 76 -3 3.6384 0.01663 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
cases_growth_daily ~ dleaving,
order = 4,
data = scc_test_big
)
## Granger causality test
##
## Model 1: cases_growth_daily ~ Lags(cases_growth_daily, 1:4) + Lags(dleaving, 1:4)
## Model 2: cases_growth_daily ~ Lags(cases_growth_daily, 1:4)
## Res.Df Df F Pr(>F)
## 1 70
## 2 74 -4 3.7127 0.008478 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ cases_growth_daily,
order = 4,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:4) + Lags(cases_growth_daily, 1:4)
## Model 2: dleaving ~ Lags(dleaving, 1:4)
## Res.Df Df F Pr(>F)
## 1 70
## 2 74 -4 2.7896 0.03281 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
cases_growth_daily ~ dleaving,
order = 5,
data = scc_test_big
)
## Granger causality test
##
## Model 1: cases_growth_daily ~ Lags(cases_growth_daily, 1:5) + Lags(dleaving, 1:5)
## Model 2: cases_growth_daily ~ Lags(cases_growth_daily, 1:5)
## Res.Df Df F Pr(>F)
## 1 67
## 2 72 -5 3.1335 0.01334 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ cases_growth_daily,
order = 5,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:5) + Lags(cases_growth_daily, 1:5)
## Model 2: dleaving ~ Lags(dleaving, 1:5)
## Res.Df Df F Pr(>F)
## 1 67
## 2 72 -5 2.8962 0.01993 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
cases_growth_daily ~ dleaving,
order = 6,
data = scc_test_big
)
## Granger causality test
##
## Model 1: cases_growth_daily ~ Lags(cases_growth_daily, 1:6) + Lags(dleaving, 1:6)
## Model 2: cases_growth_daily ~ Lags(cases_growth_daily, 1:6)
## Res.Df Df F Pr(>F)
## 1 64
## 2 70 -6 3.8148 0.002596 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ cases_growth_daily,
order = 6,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:6) + Lags(cases_growth_daily, 1:6)
## Model 2: dleaving ~ Lags(dleaving, 1:6)
## Res.Df Df F Pr(>F)
## 1 64
## 2 70 -6 2.6223 0.02457 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
cases_growth_daily ~ dleaving,
order = 7,
data = scc_test_big
)
## Granger causality test
##
## Model 1: cases_growth_daily ~ Lags(cases_growth_daily, 1:7) + Lags(dleaving, 1:7)
## Model 2: cases_growth_daily ~ Lags(cases_growth_daily, 1:7)
## Res.Df Df F Pr(>F)
## 1 61
## 2 68 -7 4.484 0.0004429 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ cases_growth_daily,
order = 7,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:7) + Lags(cases_growth_daily, 1:7)
## Model 2: dleaving ~ Lags(dleaving, 1:7)
## Res.Df Df F Pr(>F)
## 1 61
## 2 68 -7 1.7744 0.1089
grangertest(
cases_growth_daily ~ dleaving,
order = 8,
data = scc_test_big
)
## Granger causality test
##
## Model 1: cases_growth_daily ~ Lags(cases_growth_daily, 1:8) + Lags(dleaving, 1:8)
## Model 2: cases_growth_daily ~ Lags(cases_growth_daily, 1:8)
## Res.Df Df F Pr(>F)
## 1 58
## 2 66 -8 2.4973 0.02107 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ cases_growth_daily,
order = 8,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:8) + Lags(cases_growth_daily, 1:8)
## Model 2: dleaving ~ Lags(dleaving, 1:8)
## Res.Df Df F Pr(>F)
## 1 58
## 2 66 -8 2.5929 0.01706 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
cases_growth_daily ~ dleaving,
order = 9,
data = scc_test_big
)
## Granger causality test
##
## Model 1: cases_growth_daily ~ Lags(cases_growth_daily, 1:9) + Lags(dleaving, 1:9)
## Model 2: cases_growth_daily ~ Lags(cases_growth_daily, 1:9)
## Res.Df Df F Pr(>F)
## 1 55
## 2 64 -9 1.8937 0.07217 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ cases_growth_daily,
order = 9,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:9) + Lags(cases_growth_daily, 1:9)
## Model 2: dleaving ~ Lags(dleaving, 1:9)
## Res.Df Df F Pr(>F)
## 1 55
## 2 64 -9 2.2668 0.0307 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
cases_growth_daily ~ dleaving,
order = 10,
data = scc_test_big
)
## Granger causality test
##
## Model 1: cases_growth_daily ~ Lags(cases_growth_daily, 1:10) + Lags(dleaving, 1:10)
## Model 2: cases_growth_daily ~ Lags(cases_growth_daily, 1:10)
## Res.Df Df F Pr(>F)
## 1 52
## 2 62 -10 1.9322 0.06139 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
grangertest(
dleaving ~ cases_growth_daily,
order = 10,
data = scc_test_big
)
## Granger causality test
##
## Model 1: dleaving ~ Lags(dleaving, 1:10) + Lags(cases_growth_daily, 1:10)
## Model 2: dleaving ~ Lags(dleaving, 1:10)
## Res.Df Df F Pr(>F)
## 1 52
## 2 62 -10 3.5695 0.001183 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Also try ardlDlm
# 4th order lag on x only
testing_ardl4 = ardlDlm(x = scc_test_big$dleaving_mov, y = scc_test_big$d2log_cases_mov, p = 4, q = 4, remove = list(p = c(0,1,2,3), q=c()))
summary(testing_ardl4)
##
## Time series regression with "ts" data:
## Start = 5, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.0213610 -0.0014294 0.0006449 0.0016122 0.0266611
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.0009501 0.0008704 -1.092 0.279
## X.4 -0.0006286 0.0012117 -0.519 0.605
## Y.1 -0.0364700 0.1170492 -0.312 0.756
## Y.2 0.2251433 0.1197267 1.880 0.064 .
## Y.3 0.4227324 0.1001782 4.220 6.94e-05 ***
## Y.4 0.0537929 0.1099349 0.489 0.626
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.007004 on 73 degrees of freedom
## Multiple R-squared: 0.2386, Adjusted R-squared: 0.1865
## F-statistic: 4.576 on 5 and 73 DF, p-value: 0.001098
testing_ardl4 = ardlDlm(x = scc_test_big$dleaving, y = scc_test_big$d2log_cases, p = 4, q = 4, remove = list(p = c(0,1,2,3), q=c()))
summary(testing_ardl4)
##
## Time series regression with "ts" data:
## Start = 5, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.110146 -0.008604 0.001436 0.009713 0.095544
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.005328 0.003753 -1.420 0.15996
## X.4 0.001879 0.001196 1.572 0.12026
## Y.1 -0.648722 0.114527 -5.664 2.75e-07 ***
## Y.2 -0.575380 0.123212 -4.670 1.34e-05 ***
## Y.3 -0.419613 0.128582 -3.263 0.00168 **
## Y.4 -0.166485 0.093575 -1.779 0.07938 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.03261 on 73 degrees of freedom
## Multiple R-squared: 0.3428, Adjusted R-squared: 0.2978
## F-statistic: 7.616 on 5 and 73 DF, p-value: 8.376e-06
# 5th order lag on x only
testing_ardl5 = ardlDlm(x = scc_test_big$dleaving_mov, y = scc_test_big$d2log_cases_mov, p = 5, q = 5, remove = list(p = c(0,1,2,3,4), q=c()))
summary(testing_ardl5)
##
## Time series regression with "ts" data:
## Start = 6, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.0222745 -0.0009438 0.0008675 0.0020159 0.0155077
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.0014185 0.0007383 -1.921 0.0587 .
## X.5 0.0011167 0.0010129 1.102 0.2740
## Y.1 0.1180938 0.0976595 1.209 0.2306
## Y.2 -0.0262672 0.0977324 -0.269 0.7889
## Y.3 0.1833054 0.1023435 1.791 0.0775 .
## Y.4 -0.0544539 0.0932886 -0.584 0.5613
## Y.5 0.1538378 0.0919230 1.674 0.0986 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.005844 on 71 degrees of freedom
## Multiple R-squared: 0.2397, Adjusted R-squared: 0.1754
## F-statistic: 3.73 on 6 and 71 DF, p-value: 0.002788
testing_ardl5 = ardlDlm(x = scc_test_big$dleaving, y = scc_test_big$d2log_cases, p = 5, q = 5, remove = list(p = c(0,1,2,3,4), q=c()))
summary(testing_ardl5)
##
## Time series regression with "ts" data:
## Start = 6, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.101496 -0.006203 0.002880 0.011948 0.080923
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.006203 0.003331 -1.862 0.0667 .
## X.5 -0.002759 0.001062 -2.598 0.0114 *
## Y.1 -0.582852 0.101810 -5.725 2.3e-07 ***
## Y.2 -0.312538 0.119545 -2.614 0.0109 *
## Y.3 -0.311636 0.122142 -2.551 0.0129 *
## Y.4 0.199758 0.119726 1.668 0.0996 .
## Y.5 0.214992 0.083142 2.586 0.0118 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.02836 on 71 degrees of freedom
## Multiple R-squared: 0.4913, Adjusted R-squared: 0.4483
## F-statistic: 11.43 on 6 and 71 DF, p-value: 6.638e-09
# 6th order lag on x only
testing_ardl6 = ardlDlm(x = scc_test_big$dleaving_mov, y = scc_test_big$d2log_cases_mov, p = 6, q = 6, remove = list(p = c(0,1,2,3,4,5), q=c()))
summary(testing_ardl6)
##
## Time series regression with "ts" data:
## Start = 7, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.0208017 -0.0008285 0.0008797 0.0017561 0.0147256
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.0014034 0.0007658 -1.833 0.0712 .
## X.6 0.0012260 0.0010361 1.183 0.2407
## Y.1 0.1518539 0.1189929 1.276 0.2062
## Y.2 0.0062912 0.0988165 0.064 0.9494
## Y.3 0.2571029 0.0983026 2.615 0.0109 *
## Y.4 -0.0797924 0.1055060 -0.756 0.4521
## Y.5 0.1468134 0.0940073 1.562 0.1229
## Y.6 -0.1412984 0.0939433 -1.504 0.1371
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.005852 on 69 degrees of freedom
## Multiple R-squared: 0.2486, Adjusted R-squared: 0.1724
## F-statistic: 3.261 on 7 and 69 DF, p-value: 0.004746
testing_ardl6 = ardlDlm(x = scc_test_big$dleaving, y = scc_test_big$d2log_cases, p = 6, q = 6, remove = list(p = c(0,1,2,3,4,5), q=c()))
summary(testing_ardl6)
##
## Time series regression with "ts" data:
## Start = 7, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.084402 -0.003895 0.005013 0.010127 0.065618
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.0091967 0.0029874 -3.078 0.002984 **
## X.6 0.0004921 0.0010058 0.489 0.626196
## Y.1 -0.5913184 0.1037030 -5.702 2.70e-07 ***
## Y.2 -0.4276185 0.1072323 -3.988 0.000164 ***
## Y.3 -0.5602566 0.1091638 -5.132 2.52e-06 ***
## Y.4 -0.1309827 0.1114987 -1.175 0.244135
## Y.5 -0.3088103 0.1074100 -2.875 0.005367 **
## Y.6 -0.3795404 0.0769707 -4.931 5.42e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.02471 on 69 degrees of freedom
## Multiple R-squared: 0.6223, Adjusted R-squared: 0.584
## F-statistic: 16.24 on 7 and 69 DF, p-value: 1.953e-12
# 7th order lag on x only
testing_ardl7 = ardlDlm(x = scc_test_big$dleaving_mov, y = scc_test_big$d2log_cases_mov, p = 7, q = 7, remove = list(p = c(0,1,2,3,4,5,6), q=c()))
summary(testing_ardl7)
##
## Time series regression with "ts" data:
## Start = 8, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.0229043 -0.0012596 0.0003516 0.0013377 0.0111903
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.0006688 0.0006501 -1.029 0.30727
## X.7 0.0009612 0.0008592 1.119 0.26725
## Y.1 0.0922179 0.0988767 0.933 0.35435
## Y.2 0.3199778 0.0988280 3.238 0.00188 **
## Y.3 0.1504645 0.0811417 1.854 0.06809 .
## Y.4 0.1045992 0.0846600 1.236 0.22095
## Y.5 0.1618351 0.0870405 1.859 0.06738 .
## Y.6 -0.0780031 0.0785256 -0.993 0.32412
## Y.7 -0.2322415 0.0783790 -2.963 0.00421 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.004804 on 67 degrees of freedom
## Multiple R-squared: 0.4223, Adjusted R-squared: 0.3534
## F-statistic: 6.123 on 8 and 67 DF, p-value: 6.573e-06
testing_ardl7 = ardlDlm(x = scc_test_big$dleaving, y = scc_test_big$d2log_cases, p = 7, q = 7, remove = list(p = c(0,1,2,3,4,5,6), q=c()))
summary(testing_ardl7)
##
## Time series regression with "ts" data:
## Start = 8, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.057916 -0.005367 0.001760 0.007000 0.048366
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.0031241 0.0024762 -1.262 0.21146
## X.7 0.0019605 0.0007771 2.523 0.01402 *
## Y.1 -0.3011860 0.0928512 -3.244 0.00184 **
## Y.2 -0.1065801 0.0969607 -1.099 0.27561
## Y.3 -0.3201834 0.0917558 -3.490 0.00086 ***
## Y.4 0.1257686 0.0991218 1.269 0.20889
## Y.5 -0.0768311 0.0870402 -0.883 0.38055
## Y.6 -0.1397452 0.0878037 -1.592 0.11619
## Y.7 0.2323405 0.0691169 3.362 0.00128 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.01905 on 67 degrees of freedom
## Multiple R-squared: 0.7695, Adjusted R-squared: 0.742
## F-statistic: 27.96 on 8 and 67 DF, p-value: < 2.2e-16
# all
testing_ardl = ardlDlm(x = scc_test_big$dleaving_mov, y = scc_test_big$d2log_cases_mov, p = 15, q = 15)
summary(testing_ardl)
##
## Time series regression with "ts" data:
## Start = 16, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.0047279 -0.0010529 -0.0001463 0.0012977 0.0040302
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.0012875 0.0004644 -2.773 0.00875 **
## X.t 0.0017693 0.0008860 1.997 0.05343 .
## X.1 0.0009649 0.0010808 0.893 0.37791
## X.2 0.0005149 0.0009821 0.524 0.60328
## X.3 -0.0005939 0.0010123 -0.587 0.56105
## X.4 0.0002089 0.0009913 0.211 0.83429
## X.5 -0.0002293 0.0009836 -0.233 0.81701
## X.6 -0.0001107 0.0010769 -0.103 0.91872
## X.7 0.0012352 0.0011448 1.079 0.28777
## X.8 0.0001928 0.0012221 0.158 0.87555
## X.9 0.0002823 0.0010341 0.273 0.78640
## X.10 0.0020642 0.0010017 2.061 0.04662 *
## X.11 -0.0002756 0.0010028 -0.275 0.78499
## X.12 -0.0002614 0.0009913 -0.264 0.79356
## X.13 -0.0021029 0.0009363 -2.246 0.03093 *
## X.14 0.0028202 0.0009218 3.059 0.00417 **
## X.15 0.0001016 0.0009576 0.106 0.91608
## Y.1 0.1151045 0.1442918 0.798 0.43026
## Y.2 0.1030344 0.1355770 0.760 0.45222
## Y.3 0.0148629 0.1346612 0.110 0.91273
## Y.4 0.1182629 0.1258654 0.940 0.35369
## Y.5 -0.1212860 0.1159990 -1.046 0.30273
## Y.6 -0.1238151 0.1064913 -1.163 0.25261
## Y.7 -0.2529957 0.1142312 -2.215 0.03319 *
## Y.8 -0.0526962 0.0976121 -0.540 0.59262
## Y.9 0.0865780 0.0809040 1.070 0.29168
## Y.10 -0.0811094 0.0794617 -1.021 0.31419
## Y.11 -0.0081931 0.0665314 -0.123 0.90268
## Y.12 0.0450652 0.0639009 0.705 0.48520
## Y.13 0.0159947 0.0577512 0.277 0.78340
## Y.14 -0.0224993 0.0547341 -0.411 0.68346
## Y.15 -0.1010432 0.0603504 -1.674 0.10274
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.00222 on 36 degrees of freedom
## Multiple R-squared: 0.9006, Adjusted R-squared: 0.815
## F-statistic: 10.52 on 31 and 36 DF, p-value: 1.511e-10
testing_ardl = ardlDlm(x = scc_test_big$dleaving, y = scc_test_big$d2log_cases, p = 15, q = 15)
summary(testing_ardl)
##
## Time series regression with "ts" data:
## Start = 16, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.016155 -0.004704 -0.001154 0.006352 0.015003
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.293e-03 1.941e-03 -1.697 0.09841 .
## X.t 2.027e-03 6.505e-04 3.116 0.00359 **
## X.1 1.500e-03 6.940e-04 2.162 0.03736 *
## X.2 2.380e-03 6.940e-04 3.430 0.00153 **
## X.3 5.943e-04 7.676e-04 0.774 0.44383
## X.4 1.412e-03 7.068e-04 1.997 0.05338 .
## X.5 1.357e-03 7.601e-04 1.786 0.08261 .
## X.6 1.555e-04 6.938e-04 0.224 0.82397
## X.7 9.477e-04 6.627e-04 1.430 0.16132
## X.8 7.853e-05 6.578e-04 0.119 0.90563
## X.9 1.075e-03 6.735e-04 1.596 0.11919
## X.10 2.143e-03 6.861e-04 3.123 0.00352 **
## X.11 1.606e-03 7.094e-04 2.264 0.02972 *
## X.12 1.580e-03 7.038e-04 2.244 0.03104 *
## X.13 -2.323e-04 7.488e-04 -0.310 0.75823
## X.14 2.022e-03 7.047e-04 2.869 0.00684 **
## X.15 3.584e-04 7.411e-04 0.484 0.63160
## Y.1 -6.005e-01 1.615e-01 -3.718 0.00068 ***
## Y.2 -4.422e-01 1.447e-01 -3.056 0.00421 **
## Y.3 -3.223e-01 1.518e-01 -2.123 0.04070 *
## Y.4 -2.256e-01 1.451e-01 -1.554 0.12883
## Y.5 -1.919e-01 1.332e-01 -1.441 0.15830
## Y.6 -6.598e-02 1.229e-01 -0.537 0.59463
## Y.7 -7.444e-02 1.116e-01 -0.667 0.50889
## Y.8 -1.354e-01 1.057e-01 -1.281 0.20825
## Y.9 -1.528e-01 9.428e-02 -1.621 0.11380
## Y.10 -2.657e-01 9.598e-02 -2.769 0.00884 **
## Y.11 -2.394e-01 1.028e-01 -2.330 0.02554 *
## Y.12 -9.125e-02 1.042e-01 -0.876 0.38695
## Y.13 5.172e-02 9.874e-02 0.524 0.60359
## Y.14 1.025e-01 8.732e-02 1.174 0.24804
## Y.15 1.121e-01 6.369e-02 1.760 0.08683 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.0101 on 36 degrees of freedom
## Multiple R-squared: 0.7763, Adjusted R-squared: 0.5837
## F-statistic: 4.03 on 31 and 36 DF, p-value: 4.468e-05
# try with not log of cases
testing_ardl_not_log = ardlDlm(x = scc_test_big$dleaving, y = scc_test_big$d2cases, p = 15, q = 15)
summary(testing_ardl_not_log)
##
## Time series regression with "ts" data:
## Start = 16, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -14.3602 -7.9741 -0.3545 5.9212 17.7407
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.75083 1.63037 0.461 0.647909
## X.t 0.67083 0.62922 1.066 0.293465
## X.1 0.34348 0.60532 0.567 0.573945
## X.2 0.83229 0.55930 1.488 0.145437
## X.3 -0.42042 0.60086 -0.700 0.488610
## X.4 0.75119 0.59930 1.253 0.218126
## X.5 -0.62798 0.62731 -1.001 0.323479
## X.6 -0.20905 0.60994 -0.343 0.733792
## X.7 -0.80352 0.59821 -1.343 0.187608
## X.8 -0.30680 0.60241 -0.509 0.613657
## X.9 0.07638 0.58368 0.131 0.896615
## X.10 1.22742 0.57486 2.135 0.039623 *
## X.11 0.47862 0.59857 0.800 0.429184
## X.12 0.26747 0.59655 0.448 0.656584
## X.13 -0.99414 0.58608 -1.696 0.098471 .
## X.14 3.04977 0.60491 5.042 1.32e-05 ***
## X.15 0.27318 0.74302 0.368 0.715273
## Y.1 -1.31077 0.16042 -8.171 1.02e-09 ***
## Y.2 -1.47673 0.23695 -6.232 3.40e-07 ***
## Y.3 -1.76798 0.29986 -5.896 9.57e-07 ***
## Y.4 -2.03421 0.36643 -5.551 2.77e-06 ***
## Y.5 -1.86483 0.42127 -4.427 8.53e-05 ***
## Y.6 -1.71959 0.44075 -3.901 0.000402 ***
## Y.7 -1.46637 0.45636 -3.213 0.002767 **
## Y.8 -1.29391 0.45227 -2.861 0.006992 **
## Y.9 -1.12452 0.43699 -2.573 0.014337 *
## Y.10 -1.00497 0.40279 -2.495 0.017320 *
## Y.11 -0.86689 0.36373 -2.383 0.022553 *
## Y.12 -0.74504 0.32034 -2.326 0.025778 *
## Y.13 -0.56437 0.27311 -2.066 0.046038 *
## Y.14 -0.23432 0.20662 -1.134 0.264275
## Y.15 -0.11322 0.12419 -0.912 0.368015
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 11.48 on 36 degrees of freedom
## Multiple R-squared: 0.9013, Adjusted R-squared: 0.8162
## F-statistic: 10.6 on 31 and 36 DF, p-value: 1.35e-10
# pre april 15th
scc_test_pre_april15 <- scc_test_big %>%
filter(date <= "2020-04-15")
# all with this dataset
testing_ardl_pre = ardlDlm(x = scc_test_pre_april15$dleaving_mov, y = scc_test_pre_april15$d2log_cases_mov, p = 4, q = 4, remove = list(p = c(0,1,2,3)))
summary(testing_ardl_pre)
##
## Time series regression with "ts" data:
## Start = 5, End = 46
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.0193197 -0.0028480 -0.0002416 0.0022492 0.0268557
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.002229 0.001805 -1.235 0.22479
## X.4 -0.001360 0.001975 -0.689 0.49539
## Y.1 -0.044068 0.168114 -0.262 0.79471
## Y.2 0.241224 0.176943 1.363 0.18126
## Y.3 0.430497 0.143925 2.991 0.00499 **
## Y.4 0.074366 0.157262 0.473 0.63915
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.009797 on 36 degrees of freedom
## Multiple R-squared: 0.2273, Adjusted R-squared: 0.12
## F-statistic: 2.118 on 5 and 36 DF, p-value: 0.08554
testing_ardl_pre = ardlDlm(x = scc_test_pre_april15$dleaving, y = scc_test_pre_april15$d2log_cases, p = 4, q = 4, remove = list(p = c(0,1,2,3)))
summary(testing_ardl_pre)
##
## Time series regression with "ts" data:
## Start = 5, End = 46
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.101052 -0.023588 -0.003848 0.018214 0.096246
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.009740 0.007306 -1.333 0.190892
## X.4 0.002338 0.001942 1.204 0.236427
## Y.1 -0.673451 0.162609 -4.142 0.000199 ***
## Y.2 -0.614642 0.176500 -3.482 0.001322 **
## Y.3 -0.461037 0.185506 -2.485 0.017726 *
## Y.4 -0.197461 0.137284 -1.438 0.158975
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.04549 on 36 degrees of freedom
## Multiple R-squared: 0.3602, Adjusted R-squared: 0.2713
## F-statistic: 4.053 on 5 and 36 DF, p-value: 0.005082
# ardlBoundOrders
ardlBoundOrders(data = scc_test_big %>% dplyr::select(dleaving, d2log_cases), d2log_cases~dleaving)
## $p
## dleaving
## 1 13
##
## $q
## [1] 15
##
## $Stat.table
## q = 1 q = 2 q = 3 q = 4 q = 5 q = 6 q = 7
## p = 1 -311.9603 -335.3210 -330.1968 -345.8278 -366.9638 -385.0672 -376.6914
## p = 2 -320.4750 -335.3670 -332.4970 -345.0491 -367.8052 -386.9875 -378.9446
## p = 3 -330.9842 -330.9842 -337.2451 -345.8034 -367.3441 -385.7307 -377.8713
## p = 4 -322.4556 -343.6761 -343.6761 -345.1024 -365.3940 -383.7681 -375.8753
## p = 5 -321.1203 -348.4434 -346.5467 -346.5467 -363.7819 -381.8387 -374.0109
## p = 6 -347.5903 -371.9278 -369.9330 -374.0278 -374.0278 -382.2868 -374.4121
## p = 7 -340.6690 -377.1889 -376.5889 -374.8685 -375.0555 -375.0555 -373.6369
## p = 8 -354.3639 -366.6705 -364.7182 -363.3313 -374.9599 -375.8398 -375.8398
## p = 9 -350.2095 -364.4795 -362.8167 -364.8244 -372.7247 -374.0654 -379.9239
## p = 10 -369.6317 -374.5996 -372.9301 -378.0771 -379.3319 -377.5179 -377.8421
## p = 11 -345.7944 -357.9215 -356.0887 -355.1563 -364.9176 -363.8470 -365.3090
## p = 12 -358.2348 -365.2979 -365.4377 -373.5496 -373.7836 -372.2457 -371.6928
## p = 13 -392.7529 -396.8967 -402.1087 -404.1200 -404.5720 -406.0057 -405.0748
## p = 14 -388.2767 -395.8255 -398.5145 -406.6047 -404.6167 -402.7657 -400.9487
## p = 15 -401.6321 -402.1071 -410.1066 -412.1899 -411.2751 -411.8378 -411.9559
## q = 8 q = 9 q = 10 q = 11 q = 12 q = 13 q = 14
## p = 1 -384.8290 -378.0385 -373.1607 -373.3674 -376.3853 -374.1426 -382.4685
## p = 2 -386.4006 -379.9632 -373.8916 -372.9623 -375.3222 -372.7886 -380.6993
## p = 3 -387.5460 -380.4194 -375.1833 -373.8577 -378.5349 -376.0810 -381.5210
## p = 4 -387.6737 -381.6816 -376.2598 -373.6625 -377.9804 -375.6735 -383.3153
## p = 5 -387.7739 -381.9738 -375.2375 -373.2163 -376.4734 -373.7979 -383.0741
## p = 6 -386.4163 -380.7912 -374.2681 -371.6942 -375.4924 -373.1138 -384.3340
## p = 7 -384.4275 -378.8401 -372.2728 -370.3987 -373.5037 -371.2976 -382.3836
## p = 8 -385.9502 -379.2483 -373.1510 -369.4446 -371.6045 -369.5184 -380.5363
## p = 9 -379.9239 -381.0504 -375.6413 -373.6736 -373.2675 -374.9434 -383.6552
## p = 10 -377.6418 -377.6418 -377.6468 -379.1690 -376.8519 -377.2439 -389.2916
## p = 11 -369.6717 -371.1493 -371.1493 -382.4594 -381.5857 -381.1412 -397.1622
## p = 12 -371.6109 -370.0949 -381.8265 -381.8265 -379.8830 -383.2483 -399.0410
## p = 13 -403.1759 -404.1079 -406.4798 -410.7621 -410.7621 -409.4577 -410.7312
## p = 14 -399.1349 -399.6236 -401.1789 -407.4241 -405.5559 -405.5559 -409.1715
## p = 15 -409.9693 -413.1741 -414.0036 -417.9835 -418.3189 -420.3182 -420.3182
## q = 15
## p = 1 -407.8020
## p = 2 -405.8318
## p = 3 -404.0462
## p = 4 -402.3661
## p = 5 -402.8153
## p = 6 -401.2525
## p = 7 -399.9218
## p = 8 -397.9698
## p = 9 -398.4864
## p = 10 -403.0634
## p = 11 -406.7364
## p = 12 -407.4048
## p = 13 -424.1335
## p = 14 -423.2693
## p = 15 -421.5154
##
## $min.Stat
## [1] -424.1335
ARDL test might not require stationarity?
# 4th order lag on x only, original variables
testing_ardl4_2 = ardlDlm(x = scc_test_big$leaving_home_dif, y = scc_test_big$cases, p = 4, q = 4, remove = list(p = c(0,1,2,3), q=c(1,2,3,4)))
summary(testing_ardl4_2)
##
## Time series regression with "ts" data:
## Start = 5, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -7.302e-14 -6.010e-16 2.620e-16 1.336e-15 2.169e-14
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.117e-13 3.193e-15 -3.497e+01 < 2e-16 ***
## X.4 3.170e-15 1.988e-16 1.595e+01 < 2e-16 ***
## Y.1 -5.704e-16 6.813e-17 -8.373e+00 2.33e-12 ***
## Y.0 1.000e+00 6.988e-17 1.431e+16 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 9.146e-15 on 75 degrees of freedom
## Multiple R-squared: 1, Adjusted R-squared: 1
## F-statistic: 2.242e+35 on 3 and 75 DF, p-value: < 2.2e-16
# ardlBoundOrders
ardlBoundOrders(data = scc_test_big %>% dplyr::select(leaving_home_dif, cases), cases~leaving_home_dif)
## $p
## leaving_home_dif
## 1 15
##
## $q
## [1] 5
##
## $Stat.table
## q = 1 q = 2 q = 3 q = 4 q = 5 q = 6 q = 7 q = 8
## p = 1 650.2973 641.6270 635.7257 630.1846 625.1821 619.8513 607.7379 595.1162
## p = 2 641.7394 643.1781 636.6205 630.6678 625.6807 620.4806 609.3495 596.4008
## p = 3 636.3960 636.3960 638.3134 632.5898 627.6447 622.4095 611.1910 598.4001
## p = 4 629.9057 631.4029 631.4029 632.5275 626.8448 621.1701 610.8961 596.5140
## p = 5 623.7879 624.9465 626.8140 626.8140 628.8048 623.1501 612.7791 597.5042
## p = 6 615.5902 617.0908 618.9111 620.3378 620.3378 621.4988 612.8007 598.1945
## p = 7 610.4319 612.2153 613.6291 614.5274 616.3115 616.3115 614.6655 600.1745
## p = 8 604.3551 606.2545 606.9379 606.6977 607.7826 608.7765 608.7765 601.2606
## p = 9 593.0528 594.7411 595.0907 596.0349 598.0311 599.9950 590.3761 590.3761
## p = 10 587.3863 588.8949 589.3714 589.8791 591.8627 593.6131 585.9210 586.9302
## p = 11 580.3251 581.4935 582.6880 583.9051 585.6084 587.2188 579.6786 580.3983
## p = 12 574.7763 575.4343 576.9362 578.4911 580.0743 581.3619 575.0197 575.4812
## p = 13 560.2912 561.3678 559.4739 557.5347 558.9601 560.6631 560.3960 559.7173
## p = 14 534.3470 535.9078 533.5057 535.0944 529.8075 531.7172 531.1798 532.0997
## p = 15 527.2519 529.1652 525.2829 526.4698 522.1444 523.8479 523.4097 523.3561
## q = 9 q = 10 q = 11 q = 12 q = 13 q = 14 q = 15
## p = 1 588.1154 583.1995 576.7651 571.5399 566.1110 559.6037 550.7439
## p = 2 589.1846 583.7996 576.8848 571.8229 566.3383 559.7792 552.2000
## p = 3 591.1773 585.7953 578.7397 573.6795 568.2615 561.7113 554.1987
## p = 4 589.3260 584.0780 577.8451 573.0150 567.6606 560.8485 553.7741
## p = 5 589.4609 583.8559 578.1385 573.5110 568.4760 561.5495 554.1507
## p = 6 591.0250 585.6184 579.8033 575.1108 569.9874 563.0697 555.6178
## p = 7 593.0147 587.5008 581.5781 576.9497 571.8665 565.0507 557.5126
## p = 8 593.1013 587.2463 580.8001 576.2845 570.9630 564.8165 557.0980
## p = 9 592.0671 586.7181 580.8421 576.1842 570.6466 563.5534 555.0586
## p = 10 586.9302 588.3415 582.7671 578.0914 572.3785 564.7322 556.1961
## p = 11 582.3082 582.3082 584.1915 579.5716 574.1144 566.7244 558.0308
## p = 12 577.4601 579.3827 579.3827 581.3374 576.0303 568.7184 559.6443
## p = 13 561.5990 562.8479 564.6884 564.6884 564.2771 552.9216 547.1332
## p = 14 534.0976 535.3977 537.3927 538.9529 538.9529 539.5749 527.2618
## p = 15 525.2618 525.4173 527.2359 527.3845 529.2833 529.2833 529.0907
##
## $min.Stat
## [1] 522.1444
ardlBoundOrders(data = scc_test_big %>% dplyr::select(leaving_home_dif, dcases), dcases~leaving_home_dif)
## $p
## leaving_home_dif
## 1 15
##
## $q
## [1] 6
##
## $Stat.table
## q = 1 q = 2 q = 3 q = 4 q = 5 q = 6 q = 7 q = 8
## p = 1 665.0170 658.5016 651.9311 644.0689 637.0321 625.7158 612.0242 606.5196
## p = 2 658.2912 659.8454 653.6515 645.8323 638.7672 626.9381 613.5689 608.0718
## p = 3 650.2567 650.2567 651.8360 644.9846 637.5948 625.4092 613.6611 607.9474
## p = 4 645.1425 646.8280 646.8280 646.9725 639.5079 627.4010 615.5358 609.7576
## p = 5 635.8497 637.7094 639.7091 639.7091 639.9034 628.1471 616.9887 611.4031
## p = 6 621.2309 623.2291 625.1908 626.1594 626.1594 622.8506 611.5754 606.7288
## p = 7 610.4581 612.3972 614.3961 615.4649 616.7337 616.7337 611.6690 606.9245
## p = 8 618.4858 619.9157 620.8035 618.4498 615.7614 606.9768 606.9768 608.8849
## p = 9 605.2032 607.2027 608.0947 605.9398 605.3075 597.4828 595.7400 595.7400
## p = 10 599.3769 601.1624 602.6453 599.8818 598.1908 589.9821 589.3357 591.2508
## p = 11 593.8135 595.8108 597.6681 594.5465 593.2131 585.4305 584.2312 586.2021
## p = 12 583.0040 584.8137 586.3662 586.5914 586.4300 578.9500 577.8757 579.8754
## p = 13 559.2451 558.7821 560.3956 555.6145 557.1194 554.9877 555.6511 557.3983
## p = 14 536.6600 537.0920 538.6950 525.3902 527.0979 527.2312 528.4847 530.4034
## p = 15 537.3860 538.1000 539.0681 524.9590 526.7615 523.3487 524.8220 526.4190
## q = 9 q = 10 q = 11 q = 12 q = 13 q = 14 q = 15
## p = 1 600.8366 595.7055 589.1558 584.2013 578.7239 570.7727 562.7727
## p = 2 602.3376 597.3340 590.9455 585.9195 580.3933 572.5229 564.3235
## p = 3 601.8537 596.7337 591.1873 586.1220 580.1873 572.5213 563.7332
## p = 4 603.7968 598.7284 593.1797 588.0903 582.1448 574.5212 565.7221
## p = 5 605.5864 600.4345 594.5913 589.5609 583.3498 575.6893 567.6035
## p = 6 601.5609 596.3216 589.6206 584.3496 578.2406 570.8847 562.5516
## p = 7 601.9416 596.6944 590.0650 584.3657 577.7259 570.8763 562.9464
## p = 8 603.9337 598.6662 592.0493 586.3156 579.5229 572.8586 564.9062
## p = 9 597.7292 592.4713 585.5721 579.7418 571.0319 565.3131 556.3717
## p = 10 591.2508 592.9974 586.6565 580.4958 570.7486 565.6624 556.6824
## p = 11 588.1531 588.1531 588.6562 582.4936 572.5863 567.6060 558.6622
## p = 12 581.7917 583.7872 583.7872 584.2480 574.3161 569.4655 560.6430
## p = 13 559.3781 560.8696 562.2773 562.2773 558.5230 553.2529 547.6419
## p = 14 532.2670 534.0656 536.0643 537.6763 537.6763 537.5414 528.6798
## p = 15 528.3661 529.4351 531.3052 530.2919 528.3901 528.3901 529.7897
##
## $min.Stat
## [1] 523.3487
ardlBoundOrders(data = scc_test_big %>% dplyr::select(leaving_home_dif, d2cases), d2cases~leaving_home_dif)
## $p
## leaving_home_dif
## 1 14
##
## $q
## [1] 15
##
## $Stat.table
## q = 1 q = 2 q = 3 q = 4 q = 5 q = 6 q = 7 q = 8
## p = 1 708.0255 695.0573 683.2586 675.0720 638.6222 630.4959 623.7866 618.6388
## p = 2 700.9680 696.4089 684.8028 676.9067 640.6171 632.4896 625.7866 620.6376
## p = 3 687.3660 687.3660 684.3979 676.4178 641.8226 633.8255 627.3596 622.1676
## p = 4 686.2383 678.4874 678.4874 678.0445 643.2502 634.9395 628.3330 623.0245
## p = 5 680.7179 676.9110 672.2588 672.2588 645.1659 636.8969 630.3307 625.0236
## p = 6 669.9032 666.0111 664.7156 636.8153 636.8153 635.7741 629.1239 624.0158
## p = 7 657.4764 654.9395 651.9999 648.5933 630.6631 630.6631 631.1187 626.0006
## p = 8 659.8595 656.9329 653.5996 652.5340 626.0472 625.0618 625.0618 626.8475
## p = 9 652.4397 649.8297 646.0875 645.4320 619.9217 619.8871 620.3688 620.3688
## p = 10 645.4900 643.3666 639.8260 639.3664 613.5580 613.7928 613.8767 615.5973
## p = 11 637.2523 634.8880 631.7602 630.9255 605.3217 604.6739 604.8761 606.5981
## p = 12 629.9253 628.4757 625.3575 625.1048 599.8283 599.1141 598.7457 600.2866
## p = 13 612.6199 612.3519 609.7855 604.1542 587.2176 587.0744 587.4587 589.1213
## p = 14 583.0053 584.3845 568.4426 562.1544 552.3551 549.5834 546.9471 547.9067
## p = 15 574.3285 575.4507 557.0499 558.0778 545.1261 542.9300 541.0397 542.3738
## q = 9 q = 10 q = 11 q = 12 q = 13 q = 14 q = 15
## p = 1 613.0900 607.5343 602.2246 589.6840 581.7467 576.0985 569.5043
## p = 2 615.0836 609.5210 604.2234 591.6387 583.7196 578.0765 571.4873
## p = 3 616.4500 610.8603 605.7660 592.8897 584.6454 579.3115 572.5631
## p = 4 617.5113 611.9342 606.9250 594.4230 586.3950 581.1269 574.0788
## p = 5 619.4889 613.9119 608.9097 596.4163 588.3581 583.0644 576.0767
## p = 6 618.7869 613.1851 608.1755 595.9841 587.4106 581.8646 574.6438
## p = 7 620.7803 615.1715 610.1700 597.9841 589.4096 583.8613 576.5998
## p = 8 621.6175 616.0894 611.0284 598.5636 589.6877 583.9511 577.2529
## p = 9 622.1537 616.1944 611.2535 598.6515 589.8647 584.4256 577.3277
## p = 10 615.5973 617.2211 612.0947 599.4522 590.7582 585.4951 578.3009
## p = 11 608.5068 608.5068 610.5032 597.8410 589.1579 583.9348 575.2132
## p = 12 602.2501 604.2497 604.2497 598.2797 589.5414 584.0929 573.7731
## p = 13 591.1199 593.0929 588.6314 588.6314 587.5166 582.2531 573.8226
## p = 14 548.2993 549.8267 550.7611 547.5938 547.5938 547.9128 531.2568
## p = 15 542.3616 544.2568 546.0354 540.2210 539.2820 539.2820 533.1557
##
## $min.Stat
## [1] 531.2568
ardlBoundOrders(data = scc_test_big %>% dplyr::select(leaving_home_dif, d2log_cases), d2log_cases~leaving_home_dif)
## $p
## leaving_home_dif
## 1 14
##
## $q
## [1] 15
##
## $Stat.table
## q = 1 q = 2 q = 3 q = 4 q = 5 q = 6 q = 7
## p = 1 -308.6309 -330.6783 -326.9035 -342.4762 -363.8391 -385.6209 -377.3186
## p = 2 -323.3534 -338.2787 -337.0203 -345.5516 -365.7471 -384.1676 -375.9702
## p = 3 -337.0357 -337.0357 -340.0547 -345.3001 -367.3335 -386.7940 -378.4242
## p = 4 -322.4616 -347.1152 -347.1152 -346.5539 -367.4546 -386.1632 -377.7601
## p = 5 -319.8865 -344.6161 -342.6206 -342.6206 -365.7307 -384.5430 -376.0838
## p = 6 -341.4224 -364.7537 -363.6355 -368.6268 -368.6268 -382.5433 -374.0839
## p = 7 -340.8087 -378.3434 -377.5859 -375.6027 -376.7696 -376.7696 -375.0578
## p = 8 -354.5791 -366.2709 -364.3006 -362.7393 -375.7278 -375.7105 -375.7105
## p = 9 -350.4292 -363.9184 -362.3962 -364.7952 -372.1322 -371.7847 -378.4142
## p = 10 -368.8111 -373.3366 -372.4250 -376.3493 -377.4959 -375.8399 -375.5465
## p = 11 -344.8765 -356.5942 -354.8268 -354.7728 -365.0087 -363.9262 -365.3053
## p = 12 -374.0916 -377.3696 -379.7657 -384.3032 -383.1335 -381.9391 -381.5033
## p = 13 -370.0140 -370.4846 -373.4246 -379.0872 -379.0090 -379.7064 -383.9345
## p = 14 -390.6358 -399.7993 -404.4356 -412.6724 -410.7726 -408.9217 -407.4878
## p = 15 -399.5293 -399.7197 -409.3180 -413.4357 -413.2173 -412.9722 -413.1887
## q = 8 q = 9 q = 10 q = 11 q = 12 q = 13 q = 14
## p = 1 -384.1193 -378.2772 -375.7594 -373.9045 -385.5903 -387.3381 -386.0167
## p = 2 -382.8545 -376.5091 -373.9810 -372.3933 -386.4067 -391.6719 -393.2430
## p = 3 -384.4024 -378.3598 -374.2608 -371.9108 -384.8451 -390.8902 -393.3453
## p = 4 -385.5532 -378.6670 -375.5295 -372.5432 -388.0965 -393.6976 -396.4510
## p = 5 -385.9886 -379.7186 -375.7508 -372.1253 -386.3250 -392.7606 -397.2054
## p = 6 -385.7800 -380.4529 -375.2765 -372.1626 -385.7889 -392.3472 -397.2233
## p = 7 -384.4661 -379.0075 -373.7403 -370.4953 -384.0889 -390.7900 -396.5567
## p = 8 -382.4798 -377.0412 -371.8235 -369.5193 -382.1253 -388.9047 -394.6737
## p = 9 -378.4142 -377.4935 -373.0680 -368.9819 -380.9882 -387.5795 -392.6889
## p = 10 -376.4238 -376.4238 -374.7640 -372.0928 -379.5501 -388.0001 -393.1039
## p = 11 -369.8459 -370.8232 -370.8232 -377.4238 -381.6288 -387.3806 -396.6590
## p = 12 -382.0895 -380.2965 -387.4949 -387.4949 -386.6923 -391.0751 -402.1894
## p = 13 -382.8989 -380.9152 -382.8796 -382.8836 -382.8836 -394.3992 -405.4274
## p = 14 -407.1036 -407.2484 -406.7846 -411.3588 -409.4708 -409.4708 -415.8596
## p = 15 -411.1888 -413.6247 -414.2691 -418.0044 -418.2690 -420.8662 -420.8662
## q = 15
## p = 1 -406.2971
## p = 2 -407.7076
## p = 3 -405.7860
## p = 4 -404.5495
## p = 5 -402.9845
## p = 6 -403.4958
## p = 7 -401.8059
## p = 8 -400.3986
## p = 9 -398.4607
## p = 10 -398.2944
## p = 11 -402.1676
## p = 12 -405.5036
## p = 13 -406.7592
## p = 14 -422.5152
## p = 15 -421.6551
##
## $min.Stat
## [1] -422.5152
ardlBoundOrders(data = scc_test_big %>% dplyr::select(leaving_home_dif, cases_growth_daily), cases_growth_daily~leaving_home_dif)
## $p
## leaving_home_dif
## 1 14
##
## $q
## [1] 15
##
## $Stat.table
## q = 1 q = 2 q = 3 q = 4 q = 5 q = 6 q = 7 q = 8
## p = 1 455.1482 433.3015 410.9211 391.1225 383.7757 325.2643 317.0773 309.5257
## p = 2 438.6280 434.9176 407.7428 391.0769 382.9822 326.1946 318.4596 310.7772
## p = 3 429.0497 429.0497 409.2341 392.5283 384.5784 327.7080 319.1513 312.0781
## p = 4 392.5450 394.4315 394.4315 391.1500 385.6264 328.5240 320.6296 312.7734
## p = 5 405.7746 405.3126 387.3259 387.3259 384.5380 330.4206 322.5490 314.4807
## p = 6 391.2806 392.5195 370.2129 369.0508 369.0508 332.4052 324.2641 315.0921
## p = 7 345.3740 342.9651 323.9536 325.5670 326.5344 326.5344 325.7480 316.7505
## p = 8 337.8587 338.6842 331.4975 333.4971 335.2338 316.6973 316.6973 318.6885
## p = 9 340.4189 341.4923 326.0076 328.0047 329.7013 314.0485 315.4048 315.4048
## p = 10 327.3648 327.7977 313.0216 312.8989 307.4921 292.7150 294.7117 296.1535
## p = 11 308.1517 308.9899 306.0758 305.7795 307.0029 299.7709 301.6611 301.4981
## p = 12 257.2831 257.5281 258.3869 259.2928 254.8132 256.3251 255.6834 256.7494
## p = 13 245.4109 242.7546 243.4022 243.8525 236.1568 237.2273 233.9605 226.8964
## p = 14 232.6548 232.9970 223.7710 220.3923 211.3452 212.4749 214.4128 214.5589
## p = 15 229.3912 225.4392 225.1676 216.3217 209.0668 208.7627 209.3096 209.1429
## q = 9 q = 10 q = 11 q = 12 q = 13 q = 14 q = 15
## p = 1 304.6220 295.0117 289.0769 247.0309 230.1179 224.5419 205.4068
## p = 2 306.1870 296.2389 290.1364 249.0302 231.6995 225.5908 205.8936
## p = 3 306.6488 297.4570 292.0421 250.3630 233.6176 227.5879 207.7236
## p = 4 307.5799 298.2384 292.8268 251.0191 235.3207 228.6824 209.3544
## p = 5 307.5253 298.7816 294.1911 252.9023 237.3195 230.6396 211.3489
## p = 6 307.4943 297.8448 293.6651 250.4613 229.9375 226.1907 207.8137
## p = 7 309.4431 299.4357 295.6461 252.3139 230.5858 225.9718 208.4550
## p = 8 311.2740 300.7861 297.0959 253.0456 230.0256 225.7525 209.9828
## p = 9 306.5612 296.9519 294.7386 252.4602 226.8062 220.0311 209.0678
## p = 10 296.1535 297.2940 293.6219 254.4226 228.6868 221.8292 211.0547
## p = 11 296.5342 296.5342 294.1716 254.9051 230.6857 223.8239 211.0466
## p = 12 255.7004 257.0406 257.0406 254.2740 227.2089 223.0372 208.4518
## p = 13 224.5380 226.5380 225.4304 225.4304 225.9653 218.9245 206.1911
## p = 14 206.7271 207.9034 207.5801 201.3779 201.3779 203.2447 195.2065
## p = 15 210.6464 210.8963 207.1847 200.7441 195.8866 195.8866 196.5779
##
## $min.Stat
## [1] 195.2065
testing_ardl4_growth = ardlDlm(x = scc_test_big$leaving_home_dif, y = scc_test_big$cases_growth_daily, p = 5, q = 5, remove = list(p = c(0,1,2,3)))
summary(testing_ardl4_growth)
##
## Time series regression with "ts" data:
## Start = 6, End = 83
##
## Call:
## dynlm(formula = as.formula(model.text), data = data, start = 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -11.8805 -0.8568 0.0257 0.6134 12.4073
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.62782 3.16213 0.831 0.4088
## X.4 0.25736 0.13171 1.954 0.0547 .
## X.5 -0.16885 0.13355 -1.264 0.2103
## Y.1 0.25162 0.10391 2.421 0.0181 *
## Y.2 0.16684 0.10285 1.622 0.1093
## Y.3 -0.07867 0.10026 -0.785 0.4353
## Y.4 0.38987 0.08465 4.606 1.79e-05 ***
## Y.5 0.06010 0.08504 0.707 0.4820
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.344 on 70 degrees of freedom
## Multiple R-squared: 0.8532, Adjusted R-squared: 0.8385
## F-statistic: 58.12 on 7 and 70 DF, p-value: < 2.2e-16
scc_nolag <- scc_test_big %>%
dplyr::select(date, dcases, cases_growth_daily_mov, cases, log_cases, leaving_mov) %>%
cbind(`Lag` = "0")
scc_3lag <- scc_test_big %>%
dplyr::select(date, dcases, cases_growth_daily_mov, cases, log_cases, leaving3_mov) %>%
cbind(`Lag` = "3")
colnames(scc_3lag)[ncol(scc_3lag)-1] = "leaving_mov"
scc_4lag <- scc_test_big %>%
dplyr::select(date, dcases, cases_growth_daily_mov, cases, log_cases, leaving4_mov) %>%
cbind(`Lag` = "4")
colnames(scc_4lag)[ncol(scc_4lag)-1] = "leaving_mov"
scc_5lag <- scc_test_big %>%
dplyr::select(date, dcases, cases_growth_daily_mov, cases, log_cases, leaving5_mov) %>%
cbind(`Lag` = "5")
colnames(scc_5lag)[ncol(scc_5lag)-1] = "leaving_mov"
scc_6lag <- scc_test_big %>%
dplyr::select(date, dcases, cases_growth_daily_mov, cases, log_cases, leaving6_mov) %>%
cbind(`Lag` = "6")
colnames(scc_6lag)[ncol(scc_6lag)-1] = "leaving_mov"
scc_7lag <- scc_test_big %>%
dplyr::select(date, dcases, cases_growth_daily_mov, cases, log_cases, leaving7_mov) %>%
cbind(`Lag` = "7")
colnames(scc_7lag)[ncol(scc_7lag)-1] = "leaving_mov"
scc_8lag <- scc_test_big %>%
dplyr::select(date, dcases, cases_growth_daily_mov, cases, log_cases, leaving8_mov) %>%
cbind(`Lag` = "8")
colnames(scc_8lag)[ncol(scc_8lag)-1] = "leaving_mov"
scc_9lag <- scc_test_big %>%
dplyr::select(date, dcases, cases_growth_daily_mov, cases, log_cases, leaving9_mov) %>%
cbind(`Lag` = "9")
colnames(scc_9lag)[ncol(scc_9lag)-1] = "leaving_mov"
scc_10lag <- scc_test_big %>%
dplyr::select(date, dcases, cases_growth_daily_mov, cases, log_cases, leaving10_mov) %>%
cbind(`Lag` = "10")
colnames(scc_10lag)[ncol(scc_10lag)-1] = "leaving_mov"
scc_scatter <- rbind(scc_nolag, scc_3lag, scc_4lag, scc_5lag, scc_6lag, scc_7lag, scc_8lag, scc_9lag, scc_10lag)
fig_cases_growth <-
plot_ly (scc_scatter) %>%
add_trace(
x = ~leaving_mov,
y = ~cases_growth_daily_mov,
frame = ~`Lag`,
type = 'scatter',
mode = 'markers',
showlegend = F
) %>%
animation_button(visible = F) %>%
layout(xaxis = list(title = '% difference in leaving home'), yaxis = list(title = 'daily case growth rate'), title = list(title = "SCC"))
fig_cases_growth
fig_cases_abs <-
plot_ly (scc_scatter) %>%
add_trace(
x = ~leaving_mov,
y = ~dcases,
frame = ~`Lag`,
type = 'scatter',
mode = 'markers',
showlegend = F
) %>%
animation_button(visible = F) %>%
layout(xaxis = list(title = '% difference in leaving home'), yaxis = list(title = 'daily new cases'), title = list(title = "SCC"))
fig_cases_abs
For now, just plotting the metric of cumulative percent positive
smc_cases_sd_daily <- bay_sd %>%
filter(origin_census_block_group %in% smc_blockgroups$GEOID) %>%
group_by(date) %>%
summarize(total_at_home = sum(completely_home_device_count), total_devices = sum(device_count)) %>%
mutate(
percent_at_home = total_at_home*100/total_devices,
percent_leaving_home = (100 - percent_at_home),
) %>%
left_join(
smc_cumulative_cases
) %>%
filter(date >= (min(smc_cumulative_cases$date)-10))
# get the baseline percent of people at home
pre_case_growth_at_home_smc <- bay_sd %>%
filter(date < min(smc_cumulative_cases$date)) %>%
filter(origin_census_block_group %in% smc_blockgroups$GEOID) %>%
summarize(total_at_home = sum(completely_home_device_count), total_devices = sum(device_count)) %>%
mutate(percent_at_home = total_at_home*100/total_devices, percent_leaving_home = (100 - percent_at_home))
# include change in percent leaving home
smc_cases_sd_daily <- smc_cases_sd_daily %>%
mutate(leaving_home_dif = percent_leaving_home - pre_case_growth_at_home_smc$percent_leaving_home[1],
log_cases = log(cumulative_cases))
smc_test_small <- smc_cases_sd_daily %>%
mutate(
leaving_mov = movavg(leaving_home_dif, 7, type = "s"),
leaving4 = c(rep(NA,4), smc_cases_sd_daily$leaving_home_dif[1:(nrow(smc_cases_sd_daily)-4)]),
leaving4_mov = movavg(leaving4, 7, type = "s"),
leaving5 = c(rep(NA,5), smc_cases_sd_daily$leaving_home_dif[1:(nrow(smc_cases_sd_daily)-5)]),
leaving5_mov = movavg(leaving5, 7, type = "s"),
leaving6 = c(rep(NA,6), smc_cases_sd_daily$leaving_home_dif[1:(nrow(smc_cases_sd_daily)-6)]),
leaving6_mov = movavg(leaving6, 7, type = "s"),
leaving7 = c(rep(NA,7), smc_cases_sd_daily$leaving_home_dif[1:(nrow(smc_cases_sd_daily)-7)]),
leaving7_mov = movavg(leaving7, 7, type = "s"),
leaving8 = c(rep(NA,8), smc_cases_sd_daily$leaving_home_dif[1:(nrow(smc_cases_sd_daily)-8)]),
leaving8_mov = movavg(leaving8, 7, type = "s"),
leaving9 = c(rep(NA,9), smc_cases_sd_daily$leaving_home_dif[1:(nrow(smc_cases_sd_daily)-9)]),
leaving9_mov = movavg(leaving9, 7, type = "s"),
dcases = c(NA,diff(cumulative_cases)),
leaving10 = c(rep(NA,10), smc_cases_sd_daily$leaving_home_dif[1:(nrow(smc_cases_sd_daily)-10)]),
leaving10_mov = movavg(leaving10, 7, type = "s"),
dcases = c(NA,diff(cumulative_cases)),
past_cases = c(NA, smc_cases_sd_daily$cumulative_cases[1:(nrow(smc_cases_sd_daily)-1)]),
cases_growth_daily = (dcases / past_cases)*100,
cases_growth_daily_mov = movavg(cases_growth_daily, 7, type = "s")) %>%
filter(date >= "2020-03-01")
smc_test_small %>% ggplot(
aes(x = date)) +
geom_line(aes(y = leaving6, color="Leaving home")) +
geom_line(aes(y = cases_growth_daily-30, color = "Cases")) +
scale_y_continuous(sec.axis = sec_axis(~.*1+30, name = "Daily case growth rate (%)")) +
scale_colour_manual(values = c("red", "blue")) +
labs(y = "Change in percent of devices leaving home 6 days ago relative to before", x = "Date", color = "Data", title = "San Mateo County Case Growth and Social Distancing, 6 Day Lag")
# moving average
smc_test_small %>% ggplot(
aes(x = date)) +
geom_line(aes(y = leaving6_mov, color="Leaving home")) +
geom_line(aes(y = cases_growth_daily_mov-30, color = "Cases")) +
scale_y_continuous(sec.axis = sec_axis(~.*1+30, name = "Daily case growth rate (%), 7 day moving average")) +
scale_colour_manual(values = c("red", "blue")) +
labs(y = "Change in percent of devices leaving home 6 days ago relative to before, 7 day moving average", x = "Date", color = "Data", title = "San Mateo County Case Growth and Social Distancing, 6 Day Lag")
# try with percent positive cumulatively
smc_test_small %>% ggplot(
aes(x = date)) +
geom_line(aes(y = leaving6, color="Leaving home")) +
geom_line(aes(y = perc_positive_cumulative-40, color = "Cases")) +
scale_y_continuous(sec.axis = sec_axis(~.*1+40, name = "Cumulative percent of tests that are positive")) +
scale_colour_manual(values = c("red", "blue")) +
labs(y = "Change in percent of devices leaving home 6 days ago relative to before", x = "Date", color = "Data", title = "San Mateo County Positive Tests and Social Distancing, 6 Day Lag")
smc_test_small %>% ggplot(
aes(x = date)) +
geom_line(aes(y = leaving6_mov, color="Leaving home")) +
geom_line(aes(y = perc_positive_cumulative_mov-40, color = "Cases")) +
scale_y_continuous(sec.axis = sec_axis(~.*1+40, name = "Cumulative percent of tests that are positive, 7 day moving average")) +
scale_colour_manual(values = c("red", "blue")) +
labs(y = "Change in percent of devices leaving home 6 days ago relative to before, 7 day moving average", x = "Date", color = "Data", title = "San Mateo County Positive Tests and Social Distancing, 6 Day Lag")
smc_test_small %>% ggplot(
aes(x = date)) +
geom_line(aes(y = leaving6_mov, color="Leaving home")) +
geom_line(aes(y = perc_positive_movavg*100-30, color = "Cases")) +
scale_y_continuous(sec.axis = sec_axis(~.*1-30/100, name = "Percent of tests that are positive, 7 day moving average")) +
scale_colour_manual(values = c("red", "blue")) +
labs(y = "Change in percent of devices leaving home 6 days ago relative to before, 7 day moving average", x = "Date", color = "Data", title = "San Mateo County Positive Tests and Social Distancing, 6 Day Lag")
smc_nolag <- smc_test_small %>%
dplyr::select(cases_growth_daily_mov, leaving_mov) %>%
cbind(`Lag` = "0")
smc_4lag <- smc_test_small %>%
dplyr::select(cases_growth_daily_mov, leaving4_mov) %>%
cbind(`Lag` = "4")
colnames(smc_4lag)[2] = "leaving_mov"
smc_5lag <- smc_test_small %>%
dplyr::select(cases_growth_daily_mov, leaving5_mov) %>%
cbind(`Lag` = "5")
colnames(smc_5lag)[2] = "leaving_mov"
smc_6lag <- smc_test_small %>%
dplyr::select(cases_growth_daily_mov, leaving6_mov) %>%
cbind(`Lag` = "6")
colnames(smc_6lag)[2] = "leaving_mov"
smc_7lag <- smc_test_small %>%
dplyr::select(cases_growth_daily_mov, leaving7_mov) %>%
cbind(`Lag` = "7")
colnames(smc_7lag)[2] = "leaving_mov"
smc_8lag <- smc_test_small %>%
dplyr::select(cases_growth_daily_mov, leaving8_mov) %>%
cbind(`Lag` = "8")
colnames(smc_8lag)[2] = "leaving_mov"
smc_9lag <- smc_test_small %>%
dplyr::select(cases_growth_daily_mov, leaving9_mov) %>%
cbind(`Lag` = "9")
colnames(smc_9lag)[2] = "leaving_mov"
smc_10lag <- smc_test_small %>%
dplyr::select(cases_growth_daily_mov, leaving10_mov) %>%
cbind(`Lag` = "10")
colnames(smc_10lag)[2] = "leaving_mov"
smc_scatter <- rbind(smc_nolag, smc_4lag, smc_5lag, smc_6lag, smc_7lag, smc_8lag, smc_9lag, smc_10lag)
fig_cases_growth_smc <-
plot_ly (smc_scatter) %>%
add_trace(
x = ~leaving_mov,
y = ~cases_growth_daily_mov,
frame = ~`Lag`,
type = 'scatter',
mode = 'markers',
showlegend = F
) %>%
animation_button(visible = F) %>%
layout(xaxis = list(title = '% difference in leaving home'), yaxis = list(title = 'daily case growth rate'), title = list(title = "SMC"))
fig_cases_growth_smc